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Abstract 

The  goal  of  this  project  was  to  investigate  and  overcome  the  difficulties  of  extending  output- 
based  error  estimation  and  mesh  adaptation  from  steady  to  unsteady  simulations  of  aerodynam¬ 
ics  problems  relevant  to  the  U.S.  Air  Force.  These  methods  have  a  direct  impact  on  robustness, 
through  estimates  of  numerical  error  on  engineering  outputs,  and  efficiency,  through  space-time 
mesh  adaptive  methodology,  of  unsteady  simulations.  The  project  led  to  significant  advances  in 
adjoint-based  error  estimates  for  large-scale,  nonlinear,  time-dependent  problems  on  deformable 
domains,  and  in  their  use  in  driving  efficient  combined  space-time  mesh  refinement.  The  methods 
were  applied  to  several  transient  aerodynamics  problems,  culminating  with  a  three-dimensional 
adaptive  flapping-wing  simulation.  Approaches  for  reducing  the  cost  of  the  methods  were  also 
investigated,  including  a  hybridized  discontinuous  Galerkin  discretization  and  entropy-adjoint  un¬ 
steady  error  estimates.  The  findings  were  disseminated  in  a  variety  of  presentations,  conference 
proceedings,  and  archival  journal  publications. 
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1  Introduction 


Numerical  error  is  inherent  to  unsteady  Computational  Fluid  Dynamics  (CFD)  due  to  the  finite 
nature  of  computation.  To  the  extent  that  this  error  is  not  quantified,  unsteady  CFD  simulations 
are  not  robust:  using  possibly-inaccurate  results  in  a  broader  context  such  as  design  bears  risk.  As 
computations  are  performed  with  limited  resources,  of  additional  interest  is  efficiency,  which  is 
determined  by  the  cost  of  a  satisfactory  solution.  For  complex  unsteady  applications,  both  spatial 
and  temporal  resolution  are  expected  to  have  important  effects  on  numerical  error. 

Quantification  of  this  error  is  a  difficult  task  for  practitioners,  even  those  with  expert  judgment. 
Standard  error  quantification  techniques  such  as  convergence  studies  are  not  always  possible  under 
resource  constraints  and  can  depend  on  the  initial  spatial  and  temporal  resolutions.  More  robust 
are  direct  a  posteriori  estimates  of  the  output  error,  which  can  also  provide  adaptive  indicators  for 
mesh  refinement.  This  research  developed  practical  output  error  estimates  and  mesh  adaptation 
strategies  for  unsteady  flows,  with  specific  application  to  aerodynamics  problems  relevant  to  the 
Air  Force. 


2  Objectives 

The  goal  of  this  research  was  to  improve  the  robustness  and  efficiency  of  unsteady  aerodynam¬ 
ics  simulations  by  developing  an  output-based  space-time  adaptive  methodology.  The  approach 
centered  on  unsteady  discrete  adjoint-based  output  error  estimates  for  improving  robustness  and 
space-time  adaptive  mesh  refinement  for  improving  efficiency.  The  specific  objectives  were  as 
follows: 

•  To  investigate  suitable  error  measures  for  unsteady  simulations,  including  engineering  output 
error  and  surrogate  heuristics. 

•  To  extend  adjoint-based  error  estimation  theory  to  unsteady  simulations  in  a  rigorous  space- 
time  variational  setting. 

•  To  implement  an  output  error  estimation  procedure  in  a  discontinuous  Galerkin  space-time 
finite  element  discretization. 

•  To  develop  space-time  mesh  adaptation  mechanics  for  large-scale  unsteady  applications. 

•  To  demonstrate  the  proposed  error  control  framework  for  a  target  three-dimensional  appli¬ 
cation  problem. 

To  accomplish  these  objectives,  the  research  was  divided  into  three  thrusts:  high-order  discretiza¬ 
tion,  unsteady  error  estimation,  and  combined  spatial/temporal  mesh  adaptation. 


3  Accomplishments 

During  the  period  of  performance  of  this  award,  we  made  multiple  advances  in  the  proposed  re¬ 
search  objectives  related  to  unsteady  output  error  estimation  and  mesh  adaptation.  The  following 
sections  outline  the  key  contributions  with  references  to  relevant  publications. 
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3.1  High-Order  Finite  Element  Methods 

Requirements  of  high  accuracy  in  unsteady  aerodynamics  motivated  the  grounding  of  our  research 
in  high-order  spatial  and  temporal  discretizations.  The  discretization  used,  the  discontinuous 
Galerkin  (DG)  finite-element  method,  was  chosen  because  it  offers  several  key  advantages:  it 
is  unstructured,  locally  conservative,  high-order  accurate,  stable  for  convection-dominated  flows, 
low  in  dissipation  and  dispersion,  optimally  convergent,  and  amenable  to  hp  refinement.  To  make 
DG  more  robust  for  the  chosen  application  problems,  we  developed  a  constrained  pseudo-transient 
continuation  solver  [1,  5,  7]  that  enabled  us  to  simulate  challenging  cases  such  as  those  in  the 
drag-prediction  workshop  [6]  and  the  international  workshop  on  high-order  CFD  methods  [17]. 

While  most  of  our  work  was  performed  using  DG,  towards  the  end  of  the  award  we  recognized 
that  DG  methods  are  consistently  computationally  expensive:  a  discontinuous  high-order  approxi¬ 
mation  in  multiple  dimensions  creates  a  large  number  of  coupled  degrees  of  freedom  per  element. 
This  shortcoming  of  DG  motivated  the  transition  of  some  of  our  research  to  hybridized  DG  (HDG) 
methods  [9,  8]  late  in  the  award.  Figure  1  illustrates  the  primary  differences  between  DG  and 

DG  HDG 


Figure  1:  In  the  HDG  method,  additional  unknowns  on  element  interfaces  allow  elimination  of  the  element- 
interior  unknowns.  This  results  in  a  global  system  in  which  the  number  of  unknowns  scales  as  pdim~ 1  instead 
of  pdim  for  DG. 


HDG,  namely  the  introduction  of  additional  degrees  of  freedom  on  faces  that  decouple  elements 
from  each  other  and  yield  a  global  system  of  smaller  size  at  high  order  compared  to  DG.  HDG 
also  boasts  optimal  convergence  of  the  gradient  and  a  potential  local  postprocessing  technique  for 
superconvergent  solution  reconstruction. 

During  the  award  period,  we  added  new  capabilities  to  an  in-house  high-order  finite-element 
library  that  was  originally  founded  on  DG.  This  library  now  supports  both  DG  and  HDG  in  the 
same  framework.  That  is,  the  DG  and  HDG  implementations  share  functions  related  to  integra¬ 
tion,  residual  calculation  and  linearization,  and  physics;  the  only  differences  is  in  short  high-level 
functions  that  decide  which  weak-form  terms  contribute  to  the  residual.  Other  discretizations, 
such  as  the  embedded  discontinuous  Galerkin  (EDG)  wifi  now  be  easy  to  introduce  in  future  work. 
Research  efforts  supported  by  this  code  included:  robustness  improvements  for  the  Reynolds- 
averaged  Navier-Stokes  (RANS)  equations  needed  for  the  application  problems;  anisotropic  h  and 
p  refinement  for  three-dimensional  configurations;  adaptive  mesh  motion  using  an  ALE  formula¬ 
tion;  and  efficiency  and  stability  investigations  with  HDG.  Figure  2  shows  a  sample  comparison 
of  adaptive  runs  using  DG  and  HDG,  with  the  result  for  this  case  that  at  high  order,  HDG  indeed 
becomes  more  efficient  than  DG. 
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Figure  2:  .  Adaptive  simulations  of  laminar  viscous  flow  over  a  NACA  0004  airfoil  at  M  =  0.5,  Re  = 
50, 000.  A  comparison  of  DG  and  HDG  discretizations  indicates  a  benefit  of  using  HDG,  especially  at 
higher  approximation  orders. 


3.2  Error  Estimation  and  Adaptation  for  Steady-State  RANS  Flows 

As  a  precursor  to  unsteady  problems,  in  the  early  part  of  this  award  we  verified  our  steady-state 
RANS  error  estimates  and  adaptive  mesh  refinement,  and  we  extended  them  to  complex  three- 
dimensional  problems  in  preparation  for  the  unsteady  application.  The  key  added  adaptive  capabil¬ 
ity  was  Ap-rcfincmcnt  of  quadrilateral  and  hexahedral  meshes  in  two  and  three  dimensions  using  a 
discrete-choice  optimization  approach  [2,  3].  This  required  anisotropy  detection  through  element- 
patch-local  sampling  of  refinement  options  and  a  unified  output-error  driven  selection  procedure 
for  combined  p  and  anisotropic-/;  refinement  options.  We  demonstrated  the  results  on  large  three- 
dimensional  wing  and  wing-body  configurations,  using  the  RANS  equations.  A  sample  steady 
adaptive  result  is  illustrated  in  Figure  3,  which  shows  a  test  case  from  the  fifth  Drag  Prediction 
Workshop  organized  by  the  American  Institute  for  Aeronautics  and  Astronautics  in  2012  [4,  6]. 
The  case  involved  predicting  the  drag  on  an  aircraft  (wing-body)  at  transonic  cruise  flight  con¬ 
ditions.  Figure  3  shows  results  obtained  using  our  adaptive  capability,  driven  by  drag  error 
estimates.  Note  that  the  output  error  estimates  shown  in  Figure  3b  serve  as  numerical  error  bars 
that  improve  solution  robustness.  Towards  the  end  of  the  award,  we  extended  our  error  estimation 
and  adaptation  framework  to  simulations  requiring  output  trim  constraints  [16].  This  effort  also 
built  on  our  mesh  motion  capabilities  developed  as  part  of  this  award  and  discussed  in  Section  3.4. 

3.3  Output-Based  Adaptive  Space-Time  Methods 

Whereas  output-based  adaptive  methods  are  becoming  fairly  mature  for  steady  problems,  unsteady 
applications  are  less  common  because  they  pose  challenges  in  the  cost  of  error  estimates  and  in 
the  complexity  of  adaptive  mechanics.  As  part  of  this  grant,  we  extended  our  error  estimation  and 
adaptation  methods  from  steady-state  to  unsteady  simulations.  Relevant  contributions  include: 

•  Combined  temporal  and  spatial  mesh  h -refinement  with  a  static  geometry  and  mesh  [12]. 

•  Static  and  dynamic  spatial  order  refinement  on  static  geometries  [15,  10]. 

•  Combined  temporal  and  dynamic-order  spatial  refinement  on  deformable  domains  [14]. 
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(b)  Convergence  of  the  drag  output 


Figure  3:  Adaptive  results  for  the  wing  body  Common  Research  Model  (CRM),  the  subject  of  the  fifth 
AIAA  Drag  Prediction  Workshop.  Results  were  obtained  using  steady  output-based  hp  refinement  capabil¬ 
ities  of  an  in-house  (UM)  discontinuous  Galerkin  finite  element  code. 


This  work  employed  a  space-time  discontinuous  Galerkin  finite  element  discretization  using  time 
slabs  and  an  approximate  space-time  solver  [12].  For  nonlinear  problems,  the  adjoint  system  at 
each  adaptive  iteration  is  solved  backwards  in  time  using  primal  states  saved  from  the  forward  run. 
Figure  4(a)  shows  a  schematic  of  this  process.  Primal  checkpointing  was  not  needed  for  any  of 
the  runs  but  could  be  implemented  if  needed  in  the  future.  An  important  contribution  to  adaptive 
efficiency  was  in  the  form  of  a  space-time  anisotropy  measure  using  projection  of  the  line  adjoint 
to  semi-refined  spaces,  illustrated  in  Figure  4(b).  This  enabled  decisions  of  whether  spatial  or 
temporal  refinement  would  most  efficiently  reduce  the  output  error  on  a  space-time  element. 
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(a)  Unsteady  output-based  adaptive  solution  process 


(b)  Space-time  anisotropy  via  semi-refined  spaces 


Figure  4:  Schematic  of  an  adaptive  primal  and  adjoint  solution  procedure  for  output-based  unsteady  sim¬ 
ulations.  Decisions  of  refinement  in  space  and  time  are  driven  by  an  anisotropy  measure  that  projects  a 
fine-space  adjoint  down  to  semi-refined  spaces. 


Using  this  anisotropy  measure,  we  investigated  several  adaptive  schemes,  including  time  slab 
bisection,  time  node  redistribution,  static  h  spatial  mesh  refinement,  and  static/dynamic  spatial 
order  refinement.  Figure  5  shows  a  static  h  refinement  result  for  a  laminar  vortex  shedding  simula- 
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tion  [12]:  output-based  refinement  beats  uniform  and  heuristic  refinements  in  degrees  of  freedom 
for  a  strict  error  tolerance,  and  it  comes  equipped  with  an  error  estimate. 


(a)  Entropy  contours 


(b)  Output  convergence 


(c)  Output-adapted  mesh 


Figure  5:  Unsteady  output-adaptive  result  for  a  laminar  vortex  shedding  simulation  using  static  spatial 
meshes  and  DG  discretizations  in  space  and  time;  comparisons  show  performance  versus  other,  heuristic, 
indicators. 


3.4  Deformable  Domains 


In  addition  to  unsteady  problems  on  fixed  domains,  we  also  considered  domains  deforming  due  to 
motion,  such  as  that  induced  by  flapping.  These  simulations  have  far-reaching  applications,  from 
bio-inspired  flight  to  aircraft  maneuver  and  flutter  analysis.  The  runs  are  generally  computationally 
intensive  and  the  resulting  solutions  are  often  rich  in  features. 

During  the  award  period,  we  implemented  unsteady  mesh-motion  capability  into  our  solver  and 
extended  our  output-based  theory  and  methods  to  problems  involving  deformable  domains.  The 
mesh  motion  was  performed  using  an  arbitrary  Lagrangian-Eulerian  (ALE)  formulation,  illustrated 
in  Figure  6.  In  this  formulation,  error  estimation  and  adaptation  are  done  on  the  static  reference 
space.  As  such, much  of  the  output-based  theory  carries  over  in  a  straightforward  manner.  The 


(a)  Physical  space 
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(b)  Mapping 


Figure  6:  Arbitrary  Lagrangian  Eulerian  mapping  used  for  simulations  on  deforming  domains.  Error 
estimation  and  adaptation  is  done  on  the  transformed  equations  in  reference  space. 


geometric  conservation  law  (GCL)  does  require  special  attention  for  the  adjoint  solve,  error  esti¬ 
mate,  and  adaptive  indicator,  and  we  studied  the  effect  of  the  GCL  in  our  works  [13,  14].  Figure  7 
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shows  dynamic-order  adaptive  results  for  two-  and  three-dimensional  problems  with  mesh  motion. 
Output-based  refinement  clearly  beats  uniform  order  refinement  and  dynamic  order  refinement 
driven  by  a  residual  indicator.  In  our  implementation,  this  advantage  also  extends  to  computational 
time. 


(a)  2D  order  snapshot 


(c)  3D  solution  snapshots 


V)6  ‘  . 107  ‘  . 108 


Total  space-time  DOF 
(d)  3D  output  convergence 


Figure  7 :  Comparison  of  output-based  and  heuristic  space-time  refinement  strategies  for  a  two-dimensional 
problem  with  staggered  airfoils  and  a  three-dimensional  wing  undergoing  prescribed  flapping  motion  [14]. 
In  both  cases,  the  output  is  a  lift  coefficient  at/near  the  final  time  of  the  simulation. 


3.5  The  Entropy  Adjoint  Indicator 

The  entropy  adjoint  technique  for  steady-state  adaptation  was  introduced  in  collaborative  work 
prior  to  this  award.  The  key  idea  behind  this  approach  is  that  readily  available  entropy  variables 
are  in  fact  adjoint  solutions  for  an  entropy  balance  output,  and  adapting  on  an  indicator  computed 
from  this  cheap  adjoint  targets  areas  of  spurious  entropy  generation.  Work  under  the  AFOSR 
project  elucidated  connections  to  output-based  drag  adaptation  [11]  and  extended  the  technique  to 
unsteady  problems  [15]  in  a  variational  space-time  setting.  Figure  8  shows  an  adaptive  unsteady 
simulation  of  the  shallow  water  equations  using  DG  in  space  and  time:  the  entropy  adjoint  in- 
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dicator  yields  an  error  convergence  plot  most  similar  to  the  engineering-output  (force)  indicator. 
This  result  shows  that  we  can  significantly  reduce  the  cost  of  unsteady  output-based  adaptation. 
Though  the  output  is  specific  to  an  entropy  balance,  we  have  found  that  for  many  problems,  areas 
of  spurious  entropy  generation  correlate  with  areas  important  for  engineering  outputs. 


(a)  Primal  solution  (depth)  (b)  Output  convergence 


(c)  Entropy-adjoint  adapted  mesh 


Figure  8:  Adaptive  refinement  of  a  2D  shaflow-water  equation  simulation  using  an  entropy  adjoint  indica¬ 
tor:  performance  is  nearly  identical  to  an  output-based  method,  which  requires  a  separate,  backwards-in-time 
adjoint  calculation. 
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